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Abstract 

A common use of Markov Chains is tiie simulation of the seismic cycle in a fault, i.e. as a renewal 
model for the repetition of its characteristic earthquakes. This representation is consistent with 
Reid's elastic rebound theory. Here it is proved that in any one-way Markov cycle, the aperiodicity 
of the corresponding distribution of cycle lengths is always lower than one. This fact concurs with 
observations of large earthquakes in faults all over the world. 
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I. INTRODUCTION 

The elastic-rebound model is the canonical "macroscopic" theory of great earthquakes 
Q]. It states that a great earthquake will occur where large elastic strains have accumulated 
in the crust. The earthquake itself will relieve most of the strain which will then accumulate 
slowly again by a steady input of tectonic stress until the elastic strain becomes sufficiently 
large for another earthquake to ensue. The duration of this "earthquake cycle" is the ratio 
of the strain released during an earthquake to the rate of input of tectonic strain by plate 
motion. 

A corollary of the elastic-rebound model is the concept of the characteristic earthquake 
3| . Although a specific seismic fault or fault segment can produce small earthquakes as well 
as large ones, an overwhelming part of the stored elastic energy is released by the large ones, 
which tend to rupture the entire area of the fault (or fault segment). As the magnitude of an 
earthquake is related to the broken area of the fault {4], each fault (or fault segment) tends 
to produce large earthquakes of the same magnitude; and because these earthquakes release 
most of the stored elastic energy, their repetition defines the duration of the earthquake 
cycle. 

Because the Earth's crust is heterogeneous and faults are not isolated from each other 



but communicate through long-range stress-transfer mechanisms [s], the earthquake cycle 
is not periodic. So, although the elastic-rebound model is in essence deterministic, its ap- 
plication to a heterogeneous and interacting crust implies its translation into a probabilistic 
framework. Only in this way can it be used for earthquake forecasting purposes. 

Several authors have proposed probabilistic versions of the elastic-rebound model, in the 
shape of probability distribution functions (pdfs) for the duration of the earthquake cycle 
6l-lllj. The rationale of these pdfs ranges from purely statistical (e.g. Utsu [7]) to physically- 
motivated (e.g.Vazquez-Prada et al. |9|). However, due to the scarcity of registered large 
earthquakes in a specific fault (usually 4 to 10 earthquakes), the statistics upon which the 
selection of a specific pdf is based are poor. This means that different pdfs can fit the 
empirical distribution function. In any case, any pdf to be used for these purposes should 
have at least the following properties: 

1. As the time between successive earthquakes is a positive quantity, the substrate of the 
pdf must be the positive part of the real line. 



2. Very short cycles, as compared with the mean of the distribution, are extremely rare 
due to the physics of the elastic-rebound model. This means that the potential pdfs 
must be zero (or nearly so) for small values of its argument. This period of zero (or 
very small) probability at the beginning of the earthquake cycle is called stress shadow 
in seismology. 

3. The probability of cycles much longer than the mean duration is a decreasing function 
of the duration, and tends to zero as the duration tends to infinity. 

Although the complete pdf gives more information than any of its moments or other 
statistical parameters derived from it, summarizing the information content of the pdf in one 
or a few numbers is very convenient, more so when the pdf does not have a close analytical 
form or is difficult to obtain. Focusing attention on one of their statistical parameters 
benefits the comparison of the performance of different pdfs. One of these parameters is the 
coefficient of variation, a, the ratio of the standard deviation to the mean of the pdf. 

In the seismological literature the coefficient of variation is known as the aperiodicity, a 
very descriptive name when applied to the duration of the earthquake cycle: when a = the 
earthquake cycle is perfectly periodic, when < a < 1 the earthquake cycle is quasiperiodic, 
and when a > 1 the earthquake cycle is said to have a clustering of events. The case a = 1 
is particularly important because the exponential distribution has this property, and the 
exponential distribution is the pdf of an earthquake cycle where large earthquakes occur in 
time following a Poisson distribution (i.e, they are random in time). 

The predictability of a time series whose events follow a specific pdf is related to its 
aperiodicity [12]. Applied to the earthquake cycle this means that the predictability of 
the next large (characteristic) earthquake in a series is related to the aperiodicity of the pdf 
describing the duration of the cycles: aperiodicities close to zero imply greater predictability 



than aperiodicities close to one. Sykes and Menke 



13i] have calculated the aperiodicity of the 



earthquake cycle of several seismic faults. All the studied faults have aperiodicities smaller 
than 0.6, meaning that the earthquake cycle is quasiperiodic. Ellsworth et al. also 
studied the aperiodicity of the earthquake cycle in several fault segments and concluded 
that all of them are between 0.11 and 0.97. It seems that a < 1 is a property of the 
earthquake cycle in seismic faults. Can this be reproduced by simple models of single-fault 
seismicity? 
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Using the technique of Markov chains 15 1, Vazquez-Prada et ah jo] built a simple cellular 
automaton model (the minimalist model) capable of reproducing the main characteristics 
of the seismicity of individual faults. This same model was later used as a renewal model 
of seismicity 0, [l6| and applied to the forecasting of the next earthquake in the Parkfield, 
California, segment of the San Andreas Fault [11]. The aperiodicity of the earthquake cycle 
in the minimalist model is always smaller than 0.5. 

Gonzalez et al. [l3] proposed another physically-motivated Markov model of the earth- 
quake cycle (the box model). It also has characteristic-earthquake behaviour but, unlike the 
minimalist model, it only produces earthquakes of one size, that of the maximum possible 
magnitude given the area of the fault. This simplification is acceptable as large earthquakes 
release most of the stored energy. The aperiodicity of the earthquake cycle in this model is 
also lower than 1. 

A question that we want to answer in this paper is whether a general one-way Markovian 
model of the earthquake cycle can be constructed with aperiodicities larger that 1. 

The term one-way in this family of models refers to the fact that after a time step, the 
state of strain can remain either stationary or grow by a finite amount. In other words, in 
this type of model a decrease in the strain, such as could take place in a random walk model, 
is forbidden. Time increases in discrete steps and strain is also added in finite units. The 
positions of the model correspond to states of the system with progressive growing strain. 
The scheme of this model is shown in figure [H for = 6. 

The relaxation of the system through a sudden and complete loss of strain, which simu- 
lates the occurrence of an earthquake, occurs when the N^^ position of chain is reached. In 
fig. [H the relaxation is represented by the wavy line. 

This article is organized as follows: Section |IT1 contains the general form of the stochastic 
matrix of one-way Markov cycles together with the specialization to the case of the Box- 
Model and the case where all the parameters are equal. Sections IIIII and IIVI contain the 
distribution function for the Cycle Length and the two first moments of that distribution, 
respectively. The distribution function and first moments of the two particular cases men- 
tioned above are also included. Adducing that a One- Way Markov Cycle is a succession of 

independent geometric processes, a simplified form of the mean and variance are written 
in Section |Vl In Section IVIl the so called fraction of error and fraction of alarm time are 
calculated. In Section IVIIl we prove that in any one-way Markov cycle, the aperiodicity of 
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FIG. 1: Scheme of a one-way Markov cycle with = 6 



the corresponding distribution of cycle lengths is always lower than one. Finally, in Section 
IVIIII we write the conclusions. Additionally, we have considered it of interest to explicitly 
present, for a non trivial case such as = 3, how the distribution function in the case 
where all the parameters are equal tends to the Negative Binomial Distribution. This proof 
is written in the Appendix. 



II. ONE-WAY MARKOV CYCLES. TWO PARTICULAR CASES 

The Markov matrix [M] of a one-way Markov cycle of size is like this: 
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where the parameters Oj, i = 1,2,3, N, are < < 1. There are A^ positions in 
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the cycle. In each basic step, the system, with probabihty Oj, remains in the same position, 
and with probabihty (1 — Cj) jumps to the next one. This system is iUustrated in Fig. [TJ 
Denoting by [M]-^ the transpose of the Markov Matrix, the components of its eigenvector 
with eigenvahie unity, q , are: 



1 ^ 

Ci = ^ n ~ 



(1) 



TV N 



n (1 

i=l j{^i)=l 



where C is the normahzation factor. The component q is the probabihty, statistically 
speaking, of finding the system in the position i of the cycle. The value of the N components 
expressed in [1] are easily obtained from systematics. 

A particular case of this general scheme is that of the Box Model (BM) \v}\. As said in 
the introduction, this is a cellular automaton where the stochastic filling of a box represents 
the loading of elastic energy in a seismic fault. The emptying of the box after it is full is 
analogous to the generation of an earthquake. In this model, the value of the parameters 
is: 



(2) 



A^ 

Another significant particular case corresponds to the case when 

ai = a \/i (3) 

That is, all the parameters are identical. In such a Markov process, the distribution of the 
cycle length is that of a Negative Binomial Distribution (NBD) [18] where the probability 
of success is (1 — a) and successes are required. 

III. DISTRIBUTION FUNCTION FOR THE CYCLE LENGTHS 



Eliminating the element (A^, 1) of the Markov matrix in matrix [M], i. e. substituting 
(1 — Oat) by 0, we obtain a new "pruned matrix" [M]' which provides the distribution of the 
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cycle lengths, PN{n). Letter n represents the length of the cycle expressed in time steps of 
the model. Thus, n is a positive integer and PN{n) is the probability that in a system of 
size N, the cycle is completed in n steps. The formula for this probabihty is 



P^{n) = [Mr,J{l-a^) (4) 

which expresses that the probability that the length cycle is n is equal to the product of two 
factors. The first is the probability that having departed from stage 1, after n — 1 steps the 
system is in stage N. The second is the probability of passing in one step from position 
to position 1. In the first factor, the use of [M]' guarantees that in the first n—1 steps there 
has been no transition between these two positions. The computation of PN{n) provides: 



N N 



(5) 



i=i nj(^i)=i(«'i 



n = TV, + 1, . . . , oo 



This is the general form of the discrete distribution function in any one-way Markov cycle. 
This formula has been also obtained from systematics. 

As stated in the Introduction, when these systems are applied in seismicity, the fact that 
until n = N the probability of completing a cycle is null is called stress shadow. 

In the two particular cases mentioned above one obtains: 

(BM) 



(6) 



n = Ar,7V + l,...,oo 



and (NBD) 
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PtvH = (1 



n 



n 



oo 



IV. THE TWO FIRST MOMENTS 

Using eq. |5]the mean length of the cycles, //, is given by: 
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And the variance, cr^, is given by 
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V. SIMPLIFYING THE TWO FIRST MOMENTS 



In a geometric process where the probabihty of success is (1 — a), the mean and variance 
of the distribution are: 



1-a 



a 



a 



[l-af 



(14) 



Then, as a One- Way Markov Cycle is nothing more than a succession of independent 
geometric processes, the mean and variance can be written as: 

11 1 /.rX 

1 — ai 1 — 02 1 — 

and 

^ = n \^ + n \^ + ■■■ + jz ^ 16 

The reader will note the marked difference between Eqs. M and [9], which were derived 
directly from Eq. |5] and Eqs. [15] and [16] obtained by invoking the simplifying argument. 
Thus, the aperiodicity, a, is given by: 

r ai I 12 I I In 11/2 

a = ^(^-"-^ ^'-f ^'-r^'^ (17) 
— 1 h ■ ■ ■ H — 

l~ai l—a2 1— ajv 

VI. FRACTION OF ERROR AND FRACTION OF ALARM 

A simple earthquake forecasting strategy consists of turning on an alarm at a fixed value 
of n time steps after the last earthquake and maintaining the alarm until the next one. Then 



the alarm is switched off and the same strategy is repeated. 19l-l21| 

For any thinkable strategy based on the use of alarms, if an earthquake takes place when 
the alarm is on, the prediction is considered to be a success. If it takes place when the alarm 
is off, there has been a failure to predict. Then, denoting by fraction of error, /e, the number 
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of prediction failures divided by the total number of earthquakes, and the fraction of alarm, 
fa, as the ratio of the time during which the alarm is on to the total time of observation, in 
our model and using the above mentioned strategy these two functions adopt the form: 



fe{n) = = 1 - E ^(^') (18) 

n'=N n'=n+l 

Note that /e is the accumulated distribution of eq. O Performing the sum in the n' index, 
the result is: 

fe{n) = 1 - [n(l - a.)] Ef ./- ] -] (19) 



Regarding fa, 



the result is: 

oo oo 

^ifa{n) = J2 ^'Pin') -nJ2 ^(^') = (21) 

n'=n n'=n 

N N i 

f=i 1^ Uj(^i)=M - ^j) ""^y^ - ""i) (i-«i) 

VII. INEQUALITY IN THE APERIODICITY 

Knowing the value of the maximum and minimum parameters in the system, um and am, 
from eq. [T7] the aperiodicity is bounded by the following inequality: 



< < « < /^(l^^) (22) 



If all the parameters are equal to a, eq. [22] shrinks to 



or equivalently 



which is the aperiodicity of the NBD; see eqs [T2] and 
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Looking for a general inequality in a, from eq. [T7]we deduce 
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The first inequality is obvious by comparing the square of the two numerators: 

(l-ai)2 + (l-a2)2+'"^(l-a^)2- 

(26) 

. iV 1/2 1/2 

- (l-ai)2 + (l-a2)2 + ■ ■ ■ + (1 - a^)2 + 2 (1^ (1^ 
This is because all the are positive and lower than 1. 

The second inequality in eq. [25] is also obvious because each of the terms in the 
numerator is smaller than its corresponding term in the denominator. 

VIII. CONCLUSIONS 

We have calculated the form of the distribution function for the cycle length of any 
finite one-way Markov cycle. The number of independent parameters, A^, coincides with the 
number of positions in the Markov cycle. The first moments of this distribution are easily 
calculated bearing in mind that a one-way Markov cycle is nothing more than a succession 
of independent geometric processes. Thus, these moments are written as the sum of the 
mean, or variance, of the stages of the cycle. 

The above enumerated properties of the model nicely correspond to Reid's theoretical 
vision of the mechanism of how earthquakes are generated. As commented on in the Intro- 
duction, data of the time of recurrence of main shocks in faults all over the world indicate 
that the aperiodicity of those data is always lower than unity. This phenomenological fact 
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agrees with the rigorous upper bound estabhshed here: the aperiodicity of the distribution 
of cycle lengths in any one-way Markov Cycle is always lower than unity. 



APPENDIX: The Negative Binomial Distribution as a Limit: Case = 3 

In this Appendix we show explicitly that, for = 3, the limit of eq. Owhen the three 
parameters are equal is eq. [71 For simplicity in the notation, let us call CLi — (3, CL2 — 0, 
as = c. Eq. [5]for = 3 reads as follows: 



K {a — b){a — c) {b — a) (6 — c) (c — a) (c — b) 

K = (1 -a)(l -c) 

To carry out the limit, we introduce new variables x and y. 

a = xc 
b = yc 



(27) 



(28) 



The limit we seek will be implemented by tending x and y to 1. Substituting the new 
variables into eq. [271 the result is: 
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+ ^77 ^7 7T + 



c^{x — y){x—l) c^{y — x){y — 1) c^(l — — y) 

„_3 a;"-^(j/ - 1) - y^'-^jx - 1) + (x - y) 
' ^ {x-y){x-l){y-l) ^ 

Elaborating eqlSH] slightly, we obtain : 

„n-3 

-[|/(x"-^ - 1) - - 1) - y--\x - 

{x-y)[x- l)[y-l) 

Henceforth it is convenient to use the following type of polynomials: 

Pn{x) = + X"-^ + X'^-'^ + ■ ■ ■ + X + 1 

P„(l) = n + 1 
P„(x, y) = x" + x"-^i/ + x""2?/' + ■ ■ ■ + xy"-^ + 

These polynomials fulfil the so-called cyclotomic property, namely 

(x" - 1) = (x- l)P„_i(x) 
(x"-y") = (x-y)P„_i(x,i/) 

So, dividing the second factor in eq. [30] by (x — 1) we obtain 



' [2/P„_2(x)-xP„_3(x)-2/"-i] 



{x-y){y-l) 

■[2/(P„_3(x) + x"-") - xP„_3(x) - y 



„n—3 



(x-y)(2/-l) 

„n-2 „.n-2\ 



(x-?/)(|/- 1) 



x)P^_,{x)+y{x--'-y--')] 
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Now we divide the second factor of eq. [33] by {x — y) 



-^[yPn-six, y) - P„-3(x)] 



„n-3 

r[?/(a:"-^ + x^'S + x""^?/' + ■ ■ ■ + xy""" + 1/"-=^) - (x"^^ + x""" + • • ■ + x + 1)1 434) 

(2/- 1) 

„n-3 

(y- 1) 

= c"-3[x"-3 + x"-^Pi(y) + X"-5P2(2/) + ■ ■ ■ + X2p,_5(y) + xP„_4(y) + P„_3(l/)] 



Returning to eq. [271 using eq. [Ml and performing the hmit x,?/ — )■ 1, we obtain: 



hm P3(n) = (1 - cYd'-^[l + 2 + 3 + • ■ ■ + (n - 3) + (n - 2)] 



(1 - ,)a,„-.(»-lK»-2) 



(35) 



This formula coincides with eq. [7] when = 3 

P.in) = fl-c)V-3 



n — 1 
2 



(36) 



n = 3, 4, . . . , CX3 
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